Estudio Datos Vino. Depuración de variables
En este documento se exploran distintas posibilidades de estudio y depuración de un conjunto de datos con aplicación concreta al dataset de venta de vinos que se encuentra en la carpeta de Datos del material en formato cvs. Dado que la fase de “limpieza” de los datos es a crucial de cara a la obtención de modelos predictivos lo más libres posible de sesgos y problemas de estimación, es conveniente prestar atención e invertir el tiempo necesario en este proceso quizá explorando varias alternativas de las disponibles en base a la naturaleza de las distribuciones de las variables implicadas.
La fase de depuración de los datos no suele atender a “recetas” pues con cada dataset, el cuerpo nos pedirá cosas distintas. Sin embargo, trataremos de generar un proceso lo más automático posible con los aspectos generales que seguro tenemos que atacar para poder tenerlo como guía de actuación.
Estudio descriptivo de datos sobre venta de vinos
En primer lugar, importamos las dos librerías por excelencia, numpy y pandas para trabajar con conjuntos de datos en python, y leemos los datos de entrada consignando la ruta adecuada a nuestra carpeta (evidentemente hay que cambiar la ruta en cada caso, dudo que tengáis carpetar Guille/Documents…jeje)
import pandas as pd
import numpy as np
pd.set_option("display.max_rows", None, "display.max_columns", None)
# Lectura de datos
vinos = pd.read_csv('C:\\Users\\Guille\\Documents\\MineriaDatos_2022_23\\Datos\\DatosVino.csv')
vinos.head()## ID Beneficio Compra Acidez AcidoCitrico Azucar CloruroSodico \
## 0 2 515 1 0.16 -0.81 26.10 -0.425
## 1 4 585 1 2.64 -0.88 14.80 0.037
## 2 8 0 0 0.29 -0.40 21.50 0.060
## 3 11 775 1 -1.22 0.34 1.40 0.040
## 4 12 596 1 0.27 1.05 11.25 -0.007
##
## Densidad pH Sulfatos Alcohol Etiqueta CalifProductor Clasificacion \
## 0 1.02792 3.38 0.70 144.0 M 2 ***
## 1 0.99518 3.12 0.48 22.0 M 3 ***
## 2 0.99572 3.49 1.21 10.3 R 3 ?
## 3 1.03236 3.20 NaN 11.6 B 2 ***
## 4 0.99620 4.93 0.26 15.0 R 1 ?
##
## Region PrecioBotella
## 0 1.0 1.00
## 1 3.0 3.38
## 2 1.0 3.72
## 3 2.0 6.23
## 4 2.0 2.44
Este código está preparado para funcionar sobre los DatosVino. Se presenta en la siguiente tabla, la información sobre las variables contenidas en el archivo.
En la tabla de información de
variables que acompaña a los datos, podemos comprobar los tipos que
deberían tener las variables y sus posibles limitaciones sobre rangos,
categorías “oficiales” etc. Esto nos sirve de guía para ir chequeando lo
que deberíamos tener y lo que realmente encontramos (dos cosas que a
menudo no van de la mano precisamente…)
Que comprobar?
En el archivo total comprobaremos lo siguiente:
1- Tipos de variables. Todos los factores, lo son realmente en python?
2- Valores mal codificados. Hay missings no declarados tipo -1, 99999? Las categorías de las nominales son las que deben?
3- Valores fuera de rango. Hay alguna limitación sobre el rango de alguna variable que no se cumpla?
4- Variables nominales o Factores con categorías minoritarias. Las categorías con baja representación puede causar muchos problemas en los modelos por falta de base muestral para la estimación de los parámetros correspondientes a la pertenencia a esa categoría. Por ello, es conveniente echar un vistazo y recodificar las vairables uniendo categorías muy poco representativas con otras cuya unión tenga algún sentido (tienen comportamiento similar frente a la objetivo, la variable tiene caracter ordinal por lo que la unión con mayor sentido sería hacia categorías adyacentes..)
Con estas cosas ya arregladas, nos vamos a los dos grandes “caballos de batalla” de la depuración. Este proceso para la gestión de outliers y missings lo vamos a llevar a cabo sobre el archivo que llamaremos input y que contendrá solamente a los predictores sin incluir las variables objetivo continua y binaria. No es buena idea gestionar estas cosas en las objetivo puesto que las tomamos como la verdad verdadera y si tocamos algo ya estamos sesgando de alguna manera la variable de supervisión del proceso de modelización. Más adelante discutimos con mayor profundidad este hecho.
5- Outliers. Incidencia y tratamiento (pasar a missing, eliminar, winsorizar)
6- Missings. Incidencia y tratamiento (Eliminar, imputación simple por media, mediana, aleatorio, imputación por modelos)
Tipos de variables
Con la función info, podemos obtener un summary de los tipos de variables almacenados en el dataset como python lo ha entendido en la lectura.
# Información del dataset
vinos.info()## <class 'pandas.core.frame.DataFrame'>
## RangeIndex: 6365 entries, 0 to 6364
## Data columns (total 16 columns):
## # Column Non-Null Count Dtype
## --- ------ -------------- -----
## 0 ID 6365 non-null int64
## 1 Beneficio 6365 non-null int64
## 2 Compra 6365 non-null int64
## 3 Acidez 6365 non-null float64
## 4 AcidoCitrico 6365 non-null float64
## 5 Azucar 6365 non-null float64
## 6 CloruroSodico 6066 non-null float64
## 7 Densidad 6365 non-null float64
## 8 pH 6170 non-null float64
## 9 Sulfatos 5761 non-null float64
## 10 Alcohol 6365 non-null float64
## 11 Etiqueta 6365 non-null object
## 12 CalifProductor 6365 non-null int64
## 13 Clasificacion 6365 non-null object
## 14 Region 6258 non-null float64
## 15 PrecioBotella 6365 non-null float64
## dtypes: float64(10), int64(4), object(2)
## memory usage: 795.8+ KB
Se observa que no todas las variables tienen asignado el tipo correcto de datos. Identificamos factores como Compra (columna 3), Etiqueta (columna 12), Clasificación (columna 14) y Región (columna 15).
Valores únicos de las variables
Con nunique podemos obtener el número de valores distintos por cada variable para evaluar cuáles deberían ser factores.
# Número de valores distintos por variable
vinos.nunique()## ID 6365
## Beneficio 983
## Compra 2
## Acidez 659
## AcidoCitrico 523
## Azucar 1639
## CloruroSodico 1406
## Densidad 3614
## pH 436
## Sulfatos 543
## Alcohol 312
## Etiqueta 10
## CalifProductor 13
## Clasificacion 5
## Region 3
## PrecioBotella 598
## dtype: int64
En principio, claramente los factores seguros tienen menos de 10 valores distintos a excepción de Etiqueta, algo pasa ahí…Ya lo estudiaremos.
En este punto interesante charla sobre la posible dualidad continuo-categórica de la variable CalifProductor. Por una parte, tiene 13 valores, por lo que supera el umbral comentado de 10 valores distintos para ser considerada como tal. Por otra, con bajo número de valores, la linealidad con la respuesta se hace compleja, con lo que habría que comprobar la existencia de la misma. Considerando la variable como categórica, tenemos la ventaja de poder captar relaciones no lineales puesto que se modela la pertenencia a cada una de las categorías en relación a una de referencia. Sin embargo, hemos de ser conscientes del número de parámetros que consumirá en nuestro modelo \(k-1\) siendo K el número de niveles del factor.
Veamos su histograma.
import plotly.express as px
# Histograma de calificación del productor
fig = px.histogram(vinos, x="CalifProductor")
fig.show()Por tanto podemos adoptar la estrategia de mantener la variable continua por ese aspecto chi-cuadrado que no nos ha disgustado y, a su vez crear una variable categórica y posteriormente evaluar su distribución uniendo aquellas categorías minoritarias siempre con algún sentido, en este caso ordinal. Así tendremos las dos posibilidades para probar en los modelos de predicción que generemos.
Estas son posibles alternativas, en base a puro empirismo y un poco también porque se ve venir que dejar esa variable en un factor como pocas categorías que mantenga su poder predictivo es difícil tarea, decidimos mantenerla como continua.
Conversión a factor de toda variable con menos de 10 valores distintos
Utilizamos el .loc para filtrar por condición sobre valores únicos de variables y creamos la lista de columnas que deberían pasar a factor. Posteriormente, se utiliza astype para pasar a categórica.
# Lista de columnas con menos de 10 valores distintos. Potenciales factores!
to_factor = list(vinos.loc[:,vinos.nunique() < 10]);
# Podemos cambiar el tipo de todas ellas a factor de una vez
vinos[to_factor] = vinos[to_factor].astype('category')Vamos a contar el número de valores únicos de las variables numéricas por si nos hemos dejado algo por ahí.
# Número de valores distintos por variable
vinos.nunique()## ID 6365
## Beneficio 983
## Compra 2
## Acidez 659
## AcidoCitrico 523
## Azucar 1639
## CloruroSodico 1406
## Densidad 3614
## pH 436
## Sulfatos 543
## Alcohol 312
## Etiqueta 10
## CalifProductor 13
## Clasificacion 5
## Region 3
## PrecioBotella 598
## dtype: int64
En principio la única cosa rara es la variable Etiqueta que no debería tener tantos valores distintos..
Descirptivos para las variables
# Número de valores distintos por variable
vinos.describe()## ID Beneficio Acidez AcidoCitrico Azucar \
## count 6365.000000 6365.000000 6365.000000 6365.000000 6365.000000
## mean 8010.702278 452.380204 0.331214 0.314350 4718.669780
## std 4654.939139 308.380542 0.787534 0.861428 21192.546521
## min 2.000000 0.000000 -2.790000 -3.240000 -127.100000
## 25% 3980.000000 236.000000 0.130000 0.020000 0.900000
## 50% 8065.000000 480.000000 0.280000 0.310000 5.000000
## 75% 12027.000000 671.000000 0.650000 0.580000 22.600000
## max 16128.000000 1568.000000 3.680000 3.860000 99999.000000
##
## CloruroSodico Densidad pH Sulfatos Alcohol \
## count 6066.000000 6365.000000 6170.000000 5761.000000 6365.000000
## mean 0.051348 0.994204 3.202207 0.526659 16.250102
## std 0.322715 0.026417 0.678330 0.948039 25.598217
## min -1.171000 0.888090 0.540000 -3.120000 -4.500000
## 25% -0.032750 0.988245 2.960000 0.280000 9.000000
## 50% 0.046000 0.994400 3.190000 0.500000 10.500000
## 75% 0.146750 1.000600 3.460000 0.880000 12.800000
## max 1.351000 1.099240 6.050000 4.210000 150.000000
##
## CalifProductor PrecioBotella
## count 6365.000000 6365.000000
## mean 2.761508 2.610652
## std 1.319127 1.480274
## min 0.000000 1.000000
## 25% 2.000000 1.420000
## 50% 3.000000 2.190000
## 75% 3.000000 3.440000
## max 12.000000 11.440000
Distintas formas de echar el vistazo a las distribuciones de las variables, donde prestaremos atención al describe que nos informa sobre cuartiles y media, así como valores perdidos y máximos. Así, observamos que azucar tiene valores 99999 sospechosos y sulfatos 604 valores ausentes (NA), que alcohol debe tener distribución asimétrica positiva por valores posiblemente altos, de hecho es un % y no debería superar 100.
Podemos sacar los descriptivos básicos para las categóricas.
# Número de valores distintos por variable
vinos.describe(exclude=np.number)## Compra Etiqueta Clasificacion Region
## count 6365 6365 6365 6258.0
## unique 2 10 5 3.0
## top 1 R ** 3.0
## freq 4998 2380 1754 2132.0
Inspección gráfica
Podemos inspeccionar las distribuciones gráficamente para completar la exploración. En este punto, nos podemos plantear la creación de un par de funciones gráficas que nos faciliten el trabajo, de tal forma que las podamos aplicar al dataset completo para visualizarlo de un plumazo.
En primer lugar, como me gusta ver el histograma y el boxplot en relación para cada variable, vamos a definir una función que nos saque ambos gráficos jusntos. Podríamos adoptar varias formas de pensar de cara a la programación, por flexibilidad, decidimos que la función trabaje sobre una potencial columna o variable, así podremos aplicarla a todas o parte de las columnas del dataset.
import seaborn as sns
import matplotlib.pyplot as plt
def histogram_boxplot(data, xlabel = None, title = None, font_scale=2, figsize=(9,8), bins = None):
""" Boxplot and histogram combined
data: 1-d data array
xlabel: xlabel
title: title
font_scale: the scale of the font (default 2)
figsize: size of fig (default (9,8))
bins: number of bins (default None / auto)
example use: histogram_boxplot(np.random.rand(100), bins = 20, title="Fancy plot")
"""
# Definir tamaño letra
sns.set(font_scale=font_scale)
# Crear ventana para los subgráficos
f2, (ax_box2, ax_hist2) = plt.subplots(2, sharex=True, gridspec_kw={"height_ratios": (.15, .85)}, figsize=figsize)
# Crear boxplot
sns.boxplot(x=data, ax=ax_box2)
# Crear histograma
sns.histplot(x=data, ax=ax_hist2, bins=bins) if bins else sns.histplot(x=data, ax=ax_hist2)
# Pintar una línea con la media
ax_hist2.axvline(np.mean(data),color='g',linestyle='-')
# Pintar una línea con la mediana
ax_hist2.axvline(np.median(data),color='y',linestyle='--')
# Asignar título y nombre de eje si tal
if xlabel: ax_hist2.set(xlabel=xlabel)
if title: ax_box2.set(title=title, xlabel="")
# Mostrar gráfico
plt.show()
histogram_boxplot(vinos.Beneficio, bins = 20, font_scale=1, title="Distribción de Beneficio")Seguidamente vamos a implementar una función sencilla para sacar un barplot (aquí podéis indagar posibilidades de otros paquetes de visualización para customizar vuestros informes al gusto). Vamos a tirar de plotly a ver..
def cat_plot(col):
if col.dtypes == 'category':
fig = px.bar(col.value_counts())
#fig = sns.countplot(x=col)
return(fig)
# Aplicación a una variable en particular
cat_plot(vinos.Compra)cat_plot(vinos.Etiqueta)Una vez definidas las funciones para gráficos de variables numéricas y nominales, podemos integrar el proceso en una sola función que pueda ser aplicada a todas o parte de las columnas y distinga el tipo de las mismas para escoger el gráfico adecuado y mostrarlo.
def plot(col):
if col.dtypes != 'category':
print('Cont')
histogram_boxplot(col, xlabel = col.name, title = 'Distibución continua')
else:
print('Cat')
cat_plot(col)
vinos.iloc[:,1:9].apply(plot)## Cont
## Cat
## Cont
## Cont
## Cont
## Cont
## Cont
## Cont
## Beneficio None
## Compra None
## Acidez None
## AcidoCitrico None
## Azucar None
## CloruroSodico None
## Densidad None
## pH None
## dtype: object
Corrección de errores detectados
- Etiqueta tiene 10 categorías y debería tener 5. Hacemos una tabla de frecuencias de las cetegorías para indagar.
# Tabla de frecuencias
vinos.Etiqueta.value_counts()## R 2380
## M 1357
## B 1282
## r 420
## m 230
## MM 216
## b 209
## MB 191
## mm 40
## mb 40
## Name: Etiqueta, dtype: int64
Parece que hay mayúsculas y minúsculas. Habrá que unificar para que la variable tenga sentido. Lo más inmediato es pasar a mayúsculas o minúsculas todo y sencillamente se arreglará.
Podemos utilizar el apply con lambda para ver cómo qudaría el cambio (también podríamos utilizarlo para cambiarlo asignando a la variable y “pisando” la información)
# Utilizar apply con lambda para ver como quedaría el cambio
vinos['Etiqueta'].apply(lambda x: x.upper()).value_counts(sort=False)## M 1587
## R 2800
## B 1491
## MB 231
## MM 256
## Name: Etiqueta, dtype: int64
Tiene buena pinta. Podemos hacer el cambio sencillamente con operaciones sobre cadenas de caracteres de pandas .str y la función upper() para transformación en mayúsculas.
# Realizar el cambio en la propia variable y cambiar el tipo a categórica
vinos['Etiqueta'] = vinos['Etiqueta'].str.upper().astype('category')Ya que estamos con motivación, vamos a reaordenar el factor para que MM sea el primer nivel o el más bajo y MB sea el último, manteniendo la naturaleza ordinal de la variable.
# Vamos a ordenar el factor como nos gustaría tenerlo
vinos["Etiqueta"] = vinos["Etiqueta"].cat.reorder_categories(['MM','M','R','B','MB'])
# A ver como está el tema ahora
vinos.Etiqueta.value_counts()## R 2800
## M 1587
## B 1491
## MM 256
## MB 231
## Name: Etiqueta, dtype: int64
La tabla de frecuencia de value_counts() ordena por defecto las categorías por volumen de registros o frecuencia, pero siempre podemos decirle que no lo haga, con lo que mostrará el órden original del factor.
# Tabla de frecuencias sin ordenar por volumen
vinos.Etiqueta.value_counts(sort=False)## MM 256
## M 1587
## R 2800
## B 1491
## MB 231
## Name: Etiqueta, dtype: int64
- Azucar. Tiene un gráfico directamente ilegible…está claro que hay un NA no declarado con valor 99999. Hay que reemplazarlo por un verdadero missing. En python, podemos asignar un valor missing con la función nan de numpy. Al realizar el cambio con la filosofía inplace, se pisa la información original en la columna del dataset.
Cuidado con utilizar esto cuando no tenemos muy claro que el cambio es correcto!! Habría que volver atrás a datos brutos y volver a ejecutar.. que bueno, por otro lado es gratis, pero mejor para pruebas trabajar con copias o con inplace=False.
# Quitamos ese 99999 de la variable Azucar y lo pasamos a NA
vinos.Azucar.replace(99999,np.nan,inplace=True)
# Comprobamos el nuevo máximo
vinos.Azucar.max()## 141.15
- Alcohol. La tabla de info de las variables nos dice que es un % y debe etar por tanto restingido al rango 0-100. Todo valor fuera de ese rango habrá que tomarlo como valor perdido pues no podemos asignar claramente a otros valores en ppio.
# Corregimos valores fuera de rango de Alcohol (recordemos que es un %). Opción con np.where de nuumpy
#vinos.Alcohol = np.where((vinos['Alcohol'] < 0) | (vinos['Alcohol'] >100), np.nan,vinos.Alcohol)
# Opción con .loc de pandas
vinos.loc[~vinos.Alcohol.between(0, 100), "Alcohol"] = np.nan
# Comprobamos nuevos límites de la variable
vinos.Alcohol.min(),vinos.Alcohol.max()## (0.0, 26.5)
- Clasificación. Veíamos esa categoría ? que nos llama la atención. Veamos si se trata de algo residual y podemos pasar a missing para imputar posteriormente o si es una categoría con importancia en volumen y patrón frente al objetivo.
# Frecuencias de las categorías de clasificación
vinos.Clasificacion.value_counts()
# Tabla de contingencia con la variable objetivo binaria Compra## ** 1754
## ? 1680
## * 1535
## *** 1074
## **** 322
## Name: Clasificacion, dtype: int64
pd.crosstab(index=vinos['Compra'], columns=vinos['Clasificacion'])## Clasificacion * ** *** **** ?
## Compra
## 0 306 46 0 0 1015
## 1 1229 1708 1074 322 665
Reemplazar los ? por Desconocido pues representa una categoría importante en volumen y con un posible patrón interesante ante la variable objetivo. Se intuye que los vinos de clasificación desconocida tienen bastante menor probabilidad de compra!
# Reemplazar los ? por Deconocido pues representa una categoría importante en volumen y con un posible patrón interesante
vinos.Clasificacion.replace('?','Desc',inplace=True)Una vez libres de errores graves, las variables están preparadas para la gestión de outliers y missing. Para ello, es importante separa las variables objetivo y trabajar en el archivo de predictores. No es habitual tocar las variables objetivo puesto que representan nuestra verdad verdadera, son las variables de supervisión y se presuponen bien recogidas.
Imaginemos que se nos presenta la situación en la que tenemos valores missing en las objetivo, que deberíamos hacer? Pues tratar estas instancias como un conjunto de test sobre el que podríamos hacer predicciones y valorar si el modelo parece tener sentido. El problema es que no podremos evaluar la calidad de las estimaciones mediante el error cometido puesto que no tenemos su verdad verdadera.
#Indico la variableObj, el ID y las Input
# los atípicos y los missings se gestionan sólo de las input
varObjCont = vinos.Beneficio
varObjBin = vinos.Compra
imput = vinos.drop(['Beneficio','Compra'],axis=1)
imput.head()## ID Acidez AcidoCitrico Azucar CloruroSodico Densidad pH Sulfatos \
## 0 2 0.16 -0.81 26.10 -0.425 1.02792 3.38 0.70
## 1 4 2.64 -0.88 14.80 0.037 0.99518 3.12 0.48
## 2 8 0.29 -0.40 21.50 0.060 0.99572 3.49 1.21
## 3 11 -1.22 0.34 1.40 0.040 1.03236 3.20 NaN
## 4 12 0.27 1.05 11.25 -0.007 0.99620 4.93 0.26
##
## Alcohol Etiqueta CalifProductor Clasificacion Region PrecioBotella
## 0 NaN M 2 *** 1.0 1.00
## 1 22.0 M 3 *** 3.0 3.38
## 2 10.3 R 3 Desc 1.0 3.72
## 3 11.6 B 2 *** 2.0 6.23
## 4 15.0 R 1 Desc 2.0 2.44
Valores atípicos
Para facilitarnos la vida y complementar la idea que tenemos ya sobre las distribuciones de las variables, llevamos a cabo un conteo de los valores que se consideran extremos según un consenso de dos criterios distintos. En primer lugar, se distingue variable simétrica o posiblemente no, para aplicar media + 3 sd ó mediana + 8 mad, respectivamente. Recordamos en este punto que todas las medidas de dispersión basadas en la mediana o cuartiles son muy poco sensibles a la presencia de asimetría en la distribución, siendo por ello más fiables en este caso. Por otro lado, aplicamos el clásico criterio del boxplot umbrales en cuartil1 - 3IQR y cuartil3+ 3IQR.
Para aclarar, veamos un momento las asimetrías de las variables numéricas.
vinos.select_dtypes(include=np.number).apply(lambda x: x.skew())## ID 0.003339
## Beneficio 0.023872
## Acidez 0.029079
## AcidoCitrico -0.024261
## Azucar 0.022391
## CloruroSodico 0.012041
## Densidad -0.009353
## pH 0.003311
## Sulfatos -0.062373
## Alcohol 0.267099
## CalifProductor 1.663553
## PrecioBotella 1.120795
## dtype: float64
Asimetrías en valor absoluto mayores a la unidad son signo de distribución significativamente sesgada a la derecha/positiva (+) o izquierda/negativa (-).
Antes de gestionar aquellos valores detectados como outliers, valoramos la incidencia en cada variable contando el porcentaje de registros atípicos en relación al total de registros. Con ello, evaluamos si estamos ante un problema grave (tener más de un 20% de valores atípicos en alguna variable es síntoma de distribución bimodal, hay en realidad dos poblaciones jugando en la distribución) o si por el contrario, se trata de un % bajo con el que podamos lidiar sin graves consecuencias en las distribuciones de las variables.
Vamos a crear una función gestiona_outliers para facilitar el trabajo en este sentido. La idea es tener una sola función que implemente:
1- Visualización de la incidencia de outliers –> Toma de decisiones
2- Gestión de los mismos:
2.1- Winsorizar: Colapsar todo valor mayor o menor que los umbrales de cierto rango (dado por percentiles habitualmente) para igualarlo con esos umbrales. Como ya comentábamos, este método conserva la integridad del datos en el sentido de que si era muy grande por la derecha ahora estará en el extremo, pero puede provocar carga excesiva de las colas de la distribución, generando quizá algunas rarezas en picos de densidad..
2.2- Convertir en NAs: Asumiendo que los valores extremos son inaceptables (es decir, se sospecha que no puede ser tan alto o bajo por alguna razón). Podríamos pensar en que sea un fallo y por tanto decidir pasar a missing esos valores para luego gestionarlos con los valores perdidos. Esta técnica podría incurrir en errores graves de cambio del sentido de la información. Por ejemplo, en este caso un vino muy caro, pero que es un vino real y ese es su precio, es lo que hay…bueno considerándo que es erróneo, pasamos a Na e imputamos por media o valor aleatorio, quedando este valor de precio en un entorno central…tal vez hay características que generaban un patrón de precio alto que ahora serán mal aprendidas por el modelo… Hay que andarse con ojo con estas cosas..
En cualquier caso, ahí queda la función para utilizarla como se considere y tal vez incluir o desarrollar nuevas formas de tratamiento..
from scipy.stats.mstats import winsorize
# Función para gestionar outliers
def gestiona_outliers(col,clas = 'check'):
print(col.name)
# Condición de asimetría y aplicación de criterio 1 según el caso
if abs(col.skew()) < 1:
criterio1 = abs((col-col.mean())/col.std())>3
else:
criterio1 = abs((col-col.median())/col.mad())>8
# Calcular primer cuartil
q1 = col.quantile(0.25)
# Calcular tercer cuartil
q3 = col.quantile(0.75)
# Calculo de IQR
IQR=q3-q1
# Calcular criterio 2 (general para cualquier asimetría)
criterio2 = (col<(q1 - 3*IQR))|(col>(q3 + 3*IQR))
lower = col[criterio1&criterio2&(col<q1)].count()/col.count()
upper = col[criterio1&criterio2&(col>q3)].count()/col.count()
# Salida según el tipo deseado
if clas == 'check':
return(lower*100,upper*100,(lower+upper)*100)
elif clas == 'winsor':
return(winsorize(col,(lower,upper)))
elif clas == 'miss':
print('\n MissingAntes: ' + str(col.isna().sum()))
col.loc[criterio1&criterio2] = np.nan
print('MissingDespues: ' + str(col.isna().sum()) +'\n')
return(col)
# Llamada en modo check
imput.select_dtypes(include=np.number).copy().apply(lambda x: gestiona_outliers(x))## ID
## Acidez
## AcidoCitrico
## Azucar
## CloruroSodico
## Densidad
## pH
## Sulfatos
## Alcohol
## CalifProductor
## PrecioBotella
## ID Acidez AcidoCitrico Azucar CloruroSodico Densidad pH \
## 0 0.0 0.816968 0.848390 0.956307 0.807781 0.942655 0.794165
## 1 0.0 0.879811 0.785546 0.972795 0.873722 0.942655 0.696921
## 2 0.0 1.696779 1.633936 1.929101 1.681503 1.885310 1.491086
##
## Sulfatos Alcohol CalifProductor PrecioBotella
## 0 1.024128 0.000000 0.000000 0.0
## 1 0.815831 0.317885 0.094266 0.0
## 2 1.839958 0.317885 0.094266 0.0
Como pequeño resumen, lo que buscamos con el tratamiento de outliers es identificar valores verdaderamente extremos y que no sea típicos, lo cual necesariamente implica que han de ser pocos!! Y esto es muy relevante. Si son el 20% de los registros…pues tal vez estamos ante dos poblaciones distintas… habría que abordar el problema de otra forma..tal vez tratar cada población por separado.
# Crear copia para evitar pisar información
vinCont = imput.select_dtypes(include=np.number).copy()
# Aplicar la gestión de outliers en modelo winsor
vinos_winsor = vinCont.apply(lambda x: gestiona_outliers(x,clas='winsor'))
# Contemos si ha desaparecido algún resgitro o algo
# vinos_winsor.apply(lambda x: x.isna().sum()/x.count()*100)
# Veamos los valores mínimos de acidez ahora
# vinos_wins.sort_values(by='Acidez').head()
# Veamos los valores mínimos de acidez ahora
# vinos.sort_values(by='Acidez').head()## ID
## Acidez
## AcidoCitrico
## Azucar
## CloruroSodico
## Densidad
## pH
## Sulfatos
## Alcohol
## CalifProductor
## PrecioBotella
Recordemos que la búsqueda de outliers atañe en exclusiva a las variables numéricas por lo que el archivo de salida no tiene los factores. Vamos a juntar numéricas y factores en un solo archivo para valorar la incidencia de perdididos en el conjunto completo.
# Juntar columnas con join
imput_wins = vinos_winsor.join(vinos.select_dtypes(exclude=np.number))Muchas alternativas disponibles, lo importante es conocer las fortalezas y debilidades de cada uno y aplicarlo con lógica según resulte más conveniente para los datos que manejamos.
Valores perdidos
Entramos de lleno en la segunda gran gestión que debemos llevar a cabo antes de la modelización. Los archiconocidos valores perdidos. En primer lugar comentar que, llegados a este punto tenemos valores perdidos de dos fuentes distintas, por una parte los que vienen “de serie” en el dataset que podemos asociar a falta de medida en la recogida de los datos y aquí existe toda una teoría sobre los mecanismos de aparición de dichos missing, completamente al azar (MCAR), al azar(MAR) o nada de azar y existe patrón (MNAR).
Dejo por aquí un poco más de información y diversos métodos propuestos para la imputación de estos valores.
https://stefvanbuuren.name/fimd/sec-MCAR.html
Nuestro objetivo, en este limitado módulo, será valorar la incidencia de los valores perdidos y conocer los métodos usuales de imputación univariante (cada variable independientemente de las otras) con sus pros y contras.
Que podemos hacer con los missings?
- Eliminar u omitir del dataset
Eliminar “por lista” los Nas, es decir, eliminar del dataset toda observación con al menos un NA en al menos 1 variable. Se consigue liberar el dataset de perdidos (algunas técnicas y análisis no permiten su inclusión) pero se puede incurrir en los problemas que comentamos en el siguiente punto.
- Nada
Podemos obviar la presencia de valores perdidos y ya el modelo se encargará de quitarlos “por lista”, es decir, eliminar del análisis toda observación con al menos un NA en el set de variable escogido (importante este punto! si elegimos un modelo con las variables v1:v5, no se eliminarán los registros completos en estas variables pero con perdidos en v6:v10, he aquí la diferencia con la omisión de missings de primeras para toto el archivo). Esto es habitual y se puede hacer, ahora bien no está exento de peligros. Veamos.
Peligro de los missings cruzados. Imaginamos el caso en que tenemos 10 variables y 100 registros y cada una de ellas tiene un 5% de perdidos…No parece mucho con lo cual los quitamos por lista. Nuestro pensamiento es 5% pero, si se da el caso de que los NA de cada variable aparecen en registros distintos…entonces tenemos 5*10= 50% de los registros con al menos un perdido…hemos perdido la mitad de la información!!!
Sesgo por valores perdidos. El simple hecho de eliminar observaciones por el mero hecho de que presenten perdidos puede introducir un importante sesgo de selección de registros en los modelos. Imaginemos que la gente mayo tiende a no contestar ciertas preguntas en una encuesta, eliminamos y nos quedamos con muy poca gente mayor por lo que las conclusiones con toda seguridad estarán sesgadas hacia los jóvenes.
- Imputar sus valores.
No queremos exponernos a lo anterior por lo que se puede adoptar la estrategia de asignar un valor a estos datos no conocidos.
Imputación por valores centrales (media, mediana): Muy habitual asignar valores centrales ya que no alteran la distribución de las variables como tal. El gran inconveniente de este método es a subestimación de la verdadera varianza de la variable ya que estamos centrando demasiado la distribución haciendo que artificialmente su varianza se reduzca, en ocasiones muy drásticamente.
Imputación por valores aleatorios: Si no queremos centrar tanto la distribución, podemos optar por asignar al azar valores observados de las distribuciones de cada variable a los registros con NA. De esta forma, cualquier valor en el rango observado puede ser asignado a los faltantes. El gran inconveniente de esto es la dependencia del azar.
Imputación por modelos multivariantes: Muchas opciones en este apartado. Existen métodos que tienen en cuenta los valores observados de otras variables para asignar el valor más “plausible” a la variable perdida en un sentido conjunto. Una de las alternativas es generar imputaciones por un modelo de regresión por ejemplo, así para imputar Azucar utilizaremos un modelo que estime los valores de azucar en base a las demás variables (que se ajustará con los valores validos por lista) y posteriormente predeciremos los perdidos de azucar mediante este modelo generado. De esta misma forma existen modelos de imputación por random forest (missforet), vecinos más cercanos (knn), cadenas de Markov (hmisc, amelia) y gran cantidad de aproximaciones de imputación múltiple (imputar n veces y promediar). El mayor problema de estos métodos suele radicar en el posible sobreajuste a los datos de training.
El primer paso es valorar la incidencia de los valores perdidos ya que si no representan gran proporción, no existe gran peligro de cambio en la distribución de las variables con independencia del método utilizado para la imputación.
Incidencia
Vamos a actuar por dos vías, filas y columnas para contar el % de perdidos en cada registro y en cada variable. Con esta información tomaremos decisiones sobre gestión.
- Missings por variable.
#Proporción de missings por variable
imput_wins.apply(lambda x: x.isna().sum()/x.count()*100)## ID 0.000000
## Acidez 0.000000
## AcidoCitrico 0.000000
## Azucar 4.946414
## CloruroSodico 4.929113
## Densidad 0.000000
## pH 3.160454
## Sulfatos 10.484291
## Alcohol 6.491551
## CalifProductor 0.000000
## PrecioBotella 0.000000
## Compra 0.000000
## Etiqueta 0.000000
## Clasificacion 0.000000
## Region 1.709811
## dtype: float64
Controlado esto, vemos que sulfatos es la variable más peligrosa con más de un 10% de perdidos..Digamos que es la variable en la que mayores dudas tenemos respecto a la conservación de su distribución tras la imputación. Por lo demás, valores relativamente bajos.
- Missings por observación.
Calcularemos ahora el % de missings por observación y vamos a aplicar un truquiflor que a veces nos ayuda mucho a controlar las imputaciones. Se trata de generar una variable en el archivo que cuanta la proporción de perdidos que tiene cada registro. De esta forma siempre tenemos una huella de los registros con alta carga de imputaciones por si necesitamos saber algo sobre ellos en la etapa de modelado. Es más, esta será una variable que incluiremos en el modelo para valorar si puede generar un patrón de comportamiento respecto a la variable objetivo.
#Proporción de missings por observación (como una nueva columna del dataset)
imput_wins['prop_missings'] = imput_wins.apply(lambda x: x.isna().sum()/x.count()*100,axis=1)
# Valoramos distribución
imput_wins.prop_missings.describe()## count 6365.000000
## mean 2.165129
## std 3.944149
## min 0.000000
## 25% 0.000000
## 50% 0.000000
## 75% 7.142857
## max 25.000000
## Name: prop_missings, dtype: float64
Vamos a ordenar el archivo por la nueva variable creada para ver el aspecto.
imput_wins.sort_values(by='prop_missings', ascending=False).head()## ID Acidez AcidoCitrico Azucar CloruroSodico Densidad pH \
## 2781 6990 1.40 0.30 NaN NaN 1.04263 4.00
## 1344 3341 0.26 0.27 18.2 0.048 1.04227 NaN
## 5402 13598 -0.45 0.26 NaN NaN 0.99907 2.21
## 4385 11020 0.27 -0.01 NaN 0.295 0.99566 NaN
## 1352 3357 1.02 1.51 NaN NaN 0.99300 3.37
##
## Sulfatos Alcohol CalifProductor PrecioBotella Compra Etiqueta \
## 2781 NaN 15.1 3 1.30 1 R
## 1344 NaN 3.8 1 2.62 1 B
## 5402 0.68 NaN 3 3.10 1 M
## 4385 -0.46 NaN 3 2.35 1 R
## 1352 0.55 11.5 2 1.60 0 MB
##
## Clasificacion Region prop_missings
## 2781 ** 2.0 25.0
## 1344 ** NaN 25.0
## 5402 Desc 1.0 25.0
## 4385 ** 2.0 25.0
## 1352 Desc NaN 25.0
El siguiente código pretende eliminar registros y observaciones con más de la mitad de su información perdida puesto que resulta arriesgado imputar tal cantidad de datos perdidos. Es evidente que en nuestro caso no aplica y si lo ejecutamos no habrá cambio alguno.
Es importante saber que tenemos que aplicar el mismo filtro a las variables objetivo para que al unir el input depurado con ellas, los registros cuadren!
#elimino las observaciones y las variables con más de la mitad
# de datos missings (si no hay ninguna, no ejecuto este código)
input %>% filter(prop_missings< 0.5) %>% select(!(names(
prop_missingsVars)[prop_missingsVars>0.5]))-> imput
#Actualizar las observaciones de las variables objetivo
varObjBin<-varObjBin[input$prop_missings<0.5]
varObjCont<-varObjCont[input$prop_missings<0.5]Coexistencia y patrones de missings
Podemos echar un vistazo a la correlación en la existencia de missings para valorar si existe algún patrón de aparición de los mismos. En caso de que observemos patrones de aparición, podemos tirar del hilo e indagar el porqué de ese comportamiento para decidir el método más adecuado para imputar.
Hay un paquete disponible llamado missingno que implementa un par de funciones gráficas interesantes para la detección de patrones de coaparición de los valores perdidios en diferentes variables.
#conda config --add channels conda-forge
#conda install missingno
# pip install missingno
import missingno as msno
# Plot correlation heatmap of missingness
msno.matrix(imput_wins.sort_values(by='Azucar'))No se intuyen patrones de aparición de missings con facilidad.
Veamos esto en el mapa de calor.
msno.heatmap(imput_wins)Relaciones muy débiles con ese máximo de 0.1 para pH con Región. Nada destacable.
Imputaciones
Llegados a este punto, hemos comprobado que la incidencia de valores perdidos no es preocupante en general (con sulfatos con ese 10% que nos invita a estar atentos) por lo que vamos a decidirnos por la imputación de los datos de cara a mantener la mayor base muestral posible para el ajuste de los futuros modelos a los datos.
Nos puede la curiosidad y nos planteamos que pasaría si elimináramos los NAs por lista del dataset. Cuantas observaciones perderíamos? Las variables que mayor carga de perdidos presentan, son realmente relevantes para un posible modelo?
imput_wins.dropna().describe()## ID Acidez AcidoCitrico Azucar CloruroSodico \
## count 4692.000000 4692.000000 4692.000000 4692.000000 4692.000000
## mean 7996.172421 0.328282 0.324333 5.938544 0.055412
## std 4671.828755 0.766658 0.832756 34.233343 0.321277
## min 4.000000 -2.050000 -2.240000 -96.600000 -0.910000
## 25% 3917.750000 0.130000 0.030000 -1.700000 -0.029000
## 50% 8071.000000 0.280000 0.310000 4.100000 0.046000
## 75% 12026.250000 0.650000 0.580000 16.125000 0.152250
## max 16111.000000 2.680000 2.880000 140.650000 1.351000
##
## Densidad pH Sulfatos Alcohol CalifProductor \
## count 4692.000000 4692.000000 4692.000000 4692.000000 4692.000000
## mean 0.994295 3.196036 0.533186 10.586679 2.762788
## std 0.025540 0.668806 0.932729 3.498461 1.309899
## min 0.915440 1.230000 -2.150000 0.000000 0.000000
## 25% 0.988540 2.960000 0.290000 9.000000 2.000000
## 50% 0.994400 3.190000 0.500000 10.400000 3.000000
## 75% 1.000705 3.460000 0.880000 12.400000 3.000000
## max 1.072380 6.050000 4.110000 26.500000 10.000000
##
## PrecioBotella prop_missings
## count 4692.000000 4692.0
## mean 2.610217 0.0
## std 1.482378 0.0
## min 1.000000 0.0
## 25% 1.420000 0.0
## 50% 2.180000 0.0
## 75% 3.430000 0.0
## max 11.440000 0.0
Nos quedamos con 4692 registros en total para todas las variables, lo que supone un 73% de la información del archivo. Nos hemos cargado el 27% de registros..
Como veremos cuando estudiemos las relaciones con las variables objetivo, en este dataset, las variables numéricas o características químicas (todas esas distribuciones super apuntadas) no generan especial patrón con las objetivo por lo que no serán predictores relevantes. Si nos cargamos ahora información en base a sus faltantes estaremos limitando o sesgando la capacidad del modelo (que seguramente no contenga a estas variables y pueda utilizar esos registros) por disminución (además no aleatoria) de base muestral.
Conclusión, imputemos!
Vamos a explorar algunas posibilidades de los paquetes sklearn y feature_engine para imputaciones simples y multivariantes de varios tipos.
Definiremos los imputadores disponibles y posteriormente los aplicaremos a los datos.
import sklearn.impute as skl_imp
from sklearn.experimental import enable_iterative_imputer
# Moda: Solo nominales
imputer_moda = skl_imp.SimpleImputer(strategy='most_frequent', missing_values=np.nan)
# knn: Solo numéricas
imputer_knn = skl_imp.KNNImputer(n_neighbors=3)
# Chain equations: solo numéricas
imputer_itImp = skl_imp.IterativeImputer(max_iter=10, random_state=0)
import feature_engine.imputation as fe_imp
# Aleatoria: numéricas y nominales
imputer_rand = fe_imp.RandomSampleImputer()
# Mediana: solo numéricas
imputer_median = fe_imp.MeanMedianImputer(imputation_method='median')
# Media: solo nominales
imputer_mean = fe_imp.MeanMedianImputer(imputation_method='mean')Separaremos el dataset de nuevo en continuas y categóricas para aplicar los métodos que correspondan.
imput_wins_cont = imput_wins.select_dtypes(include=np.number)
imput_wins_cat = imput_wins.select_dtypes(exclude=np.number)Posibilidades para las numéricas:
1- Nivel univariante
Comenzamos proponiendo métodos de imputación univariante que se basan en la información exclusiva de la distribución de la propia variable sin mirar a sus compañeras. En este plano, opciones ya comentadas como media, mediana o aleatorio.
# Media
vinos_winsor_mean_imputed = imputer_mean.fit(imput_wins_cont).transform(imput_wins_cont)
# Mediana
vinos_winsor_median_imputed = imputer_median.fit(imput_wins_cont).transform(imput_wins_cont)2- Nivel multivariante
Podemos explorar un par de posibilidades que se utilizan en el mundillo. Por una parte, las imputaciones basadas en la técnica de los k vecinos más cercanos (digamos que es un método espacial) que asocia valores promedio de los k vecinos más cercanos a la observación perdida pero jugando en el espación R^k siendo k el número de variables. Es decir, identifica registros parecidos en general en todas las variables (que estén cerca en ese hiperespacio), lo cual es una buena idea. El mayor problema que tiene es la dependencia del valor k (no hay consenso en este aspecto y depende de los datos) y el posible sobreajuste a los datos de training como cualquier otro método multi.
# Fit/transform
imput_wins_knn_imputed = pd.DataFrame(imputer_knn.fit_transform(imput_wins_cont),columns=imput_wins_cont.columns)
imput_wins_itImp_imputed = pd.DataFrame(imputer_itImp.fit_transform(imput_wins_cont),columns=imput_wins_cont.columns)Posibilidades para las nominales:
Para las variables categóricas podemos utilizar moda (la categoría más representada) o aleatorio.
imput_wins_moda_imputed = pd.DataFrame(imputer_moda.fit_transform(imput_wins_cat),columns=imput_wins_cat.columns)Aplicar random a nivel general ya que acepta todo tipo.
imput_wins_rand_imputed = imputer_rand.fit(imput_wins).transform(imput_wins)
imput_wins_rand_imputed.describe()## ID Acidez AcidoCitrico Azucar CloruroSodico \
## count 6365.000000 6365.000000 6365.000000 6365.000000 6365.000000
## mean 8010.702278 0.330902 0.315136 5.801131 0.051932
## std 4654.939139 0.768507 0.841088 33.980154 0.318607
## min 2.000000 -2.050000 -2.240000 -96.600000 -0.910000
## 25% 3980.000000 0.130000 0.020000 -1.600000 -0.032000
## 50% 8065.000000 0.280000 0.310000 4.000000 0.046000
## 75% 12027.000000 0.650000 0.580000 15.900000 0.147000
## max 16128.000000 2.680000 2.880000 141.150000 1.351000
##
## Densidad pH Sulfatos Alcohol CalifProductor \
## count 6365.000000 6365.000000 6365.000000 6365.000000 6365.000000
## mean 0.994207 3.204671 0.542364 10.625656 2.760094
## std 0.025614 0.671902 0.932792 3.547211 1.310444
## min 0.915440 1.230000 -2.150000 0.000000 0.000000
## 25% 0.988245 2.960000 0.290000 9.000000 2.000000
## 50% 0.994400 3.190000 0.500000 10.400000 3.000000
## 75% 1.000600 3.470000 0.890000 12.400000 3.000000
## max 1.072380 6.050000 4.210000 26.500000 10.000000
##
## PrecioBotella prop_missings
## count 6365.000000 6365.000000
## mean 2.610652 2.165129
## std 1.480274 3.944149
## min 1.000000 0.000000
## 25% 1.420000 0.000000
## 50% 2.190000 0.000000
## 75% 3.440000 7.142857
## max 11.440000 25.000000
Toma de decisiones y guardado del archivo depurado
# Agregar variables objetivo al input ya limpio
vinos_wins_rand_imputed = pd.concat([imput_wins_rand_imputed, varObjCont,varObjBin], axis=1)
# Guardar archivo
vinos_wins_rand_imputed.to_csv('C:\\Users\\Guille\\Documents\\MineriaDatos_2022_23\\PARTE I_Depuracion y Regresiones\\Dia1_MDDepuracion\\DatosVinoDep_winsRand.csv')Ya tenemos los datos depurados para poder empezar con el modelado. Es importante saber que a la hora de modelar utilizaremos este nuevo conjunto datosVinoDep_metodos y no el original, que para eso nos lo hemos trabajado.
Notas:
Debido a la elección del método de imputación por valores aleatorios, cada ejecución de este script generará un conjunto de datos depurados con ligeras diferencia en distribución de las variables sometidas a dicho proceso. Esto es natural y, dado que la verdadera distribución de los valores de los perdidos se deconoce en teoría, cualquiera de las soluciones sería plausible. Este comportamiento provoca que el posterior proceso de modelización con base en los datos depurados presente valores dependientes de la ejecución concreta de este primer paso. La idea es que nuestra limpieza de los datos sea relativamente robusta, por lo que los cambios deberían ser muy ligeros y no producirse cosas como cambio radical de influencia de las variables en los modelos o muy distintas capacidades de predicción sobre los datos de test. Si esto pasa, sospecharemos que nuestro proceso preliminar produce sesgos peligrosos y tendremos que repasar en especial el tratamiento de outliers y su conversión a valores perdidos (son lo que más influencia tienen en este hecho!)
Es evidente que muchas tipos de tratamiento de variables son posibles y normalmente es buena política “jugar” lo máximo con ellos y proponer distintas recodificaciones para aumentar nuestras posibilidades en la búsqueda de esos patrones ocultos que relacionan la variabilidad de las variables y que serán descubiertos por el modelo para ofrecer predicciones decentes.
Antes de lanzarse a las imputaciones es conveniente pensar si existen relaciones entre variables de manera natural, por ejemplo cuando dos continuas tienen una relación lineal fortísima digamos 0.9 de coeficiente de correlación r de Pearson, podría ser muy buena idea asignar directamente valores faltantes de una de ellas cuando la otra presenta valor conocido, así podemos reducir la incertidumbre y el sesgo de las imputaciones. De igual forma, relaciones transitivas tipo A ~ B*C (el producto de B y C presenta una cotrrelación alta con la variable A). Oye pues cuando A está perdida y ese registro tiene los valores de B y C conocido, en lugar de imputarle su valor de A, asignamos directamente.
En cuanto a código y ejecuciones en RMD, por las dudas. Los chunks contienen el código y pueden ser ejecutados por separado en su totalidad con el boton del play (parte superior derecha del chunk), línea a línea con ctrl+enter o el botón Run. Hemos de saber que a la hora de compilar el documento, las rutas deben ser válidas y contener los elementos que se requieren! El tipo de documento que se genera se controla en la cabecera y hay muchísimas posibilidades tanto en html como en word o pdf. Os animo a investigar formatos de vuestro agrado.